Staphylococcus aureus isolates from children with clinically differentiated osteomyelitis exhibit distinct transcriptomic signatures

There is substantial genomic heterogeneity among Staphylococcus aureus isolates of children with acute hematogenous osteomyelitis (AHO) but transcriptional behavior of clinically differentiated strains has not been previously described. This study evaluates transcriptional activity of S. aureus isolates of children with AHO that may regulate metabolism, biosynthesis, or virulence during bacterial growth and pathogenesis. In vitro growth kinetics were compared between three S. aureus clinical isolates from children with AHO who had mild, moderate, and severe illness. Total RNA sequencing was performed for each isolate at six separate time points throughout the logarithmic phase of growth. The NASA RNA-Sequencing Consensus Pipeline was used to identify differentially expressed genes allowing for 54 comparisons between the three isolates during growth. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment pathways were used to evaluate transcriptional variation in metabolism, biosynthesis pathways and virulence potential of the isolates. The S. aureus isolates demonstrated differing growth kinetics under standardized conditions with the mild isolate having higher optical densities with earlier and higher peak rates of growth than that of the other isolates (p<0.001). Enrichment pathway analysis established distinct transcriptional signatures according to both sampling time and clinical severity. Moderate and severe isolates demonstrated pathways of bacterial invasion, S. aureus infection, quorum sensing and two component systems. In comparison, the mild strain favored biosynthesis and metabolism. These findings suggest that transcriptional regulation during the growth of S. aureus may impact the pathogenetic mechanisms involved in the progression of severity of illness in childhood osteomyelitis. The clinical isolates studied demonstrated a tradeoff between growth and virulence. Further investigation is needed to evaluate these transcriptional pathways in an animal model or during active clinical infections of children with AHO.


Introduction
Children with acute hematogenous osteomyelitis (AHO) caused by Staphylococcus aureus demonstrate a wide range of illnesses, ranging from mild to severe [1,2]. Some children have short hospitalizations without surgery, while others develop septic shock and have prolonged bacteremia, requiring intensive care and surgical source control during prolonged hospitalizations [1,2]. The underlying pathogenetic mechanisms which lead to this wide spectrum of clinical phenotypes of AHO have yet to be established with respect to either the host or pathogen. Previous studies have evaluated the host gene expression from blood samples of children with S. aureus AHO and observed up-regulation of innate immunity (mainly neutrophil activity) and down regulation of adaptive immunity (T cells, B cells, and NK cells) during the period of acute infection [3][4][5]. However, it is unclear how these activities potentiate differentiation toward mild or severe illness. From the perspective of the pathogen, a high-resolution transcriptomic analysis of S. aureus in a mouse hematogenous model identified differentially expressed genes in acute and chronic osteomyelitis [6]. Genes that mediate metabolic adaptation, immune evasion, and replication appear to drive acute osteomyelitis. Whereas in chronic osteomyelitis the pathogen switches its transcriptional response to a persistence mode driven by nutritional deficiencies [6]. This suggests that transcriptional differences of the pathogen may play an important role in the development and progression of AHO. This theory is further supported by substantial genomic heterogeneity which has been observed among S. aureus isolates obtained from children with AHO within a single community [7]. A phylogenetic analysis of seventy-one clinical isolates from children with a wide spectrum of clinical phenotypes, ranging from mild to severe, demonstrated correlation of the genetic distances between isolates and the severity of illness scores of the affected children [8]. This study analyzes and demonstrates important differences in the growth kinetics and transcriptional signatures of S. aureus isolates obtained from children with clinically differentiated AHO.

Bacterial strain selection
This study was conducted following Institutional Review Board (IRB) approval: IRB Number STU 032013-022. The S. aureus isolates were identified from a repository of clinical strains previously evaluated with whole genome sequencing [7,8]. Organisms were selected from opposite ends of the phylogenetic map to ensure maximum genetic distance between isolates and wide clinical differentiation based on the severity of illness scores of the affected children. Additionally, one pathogen was chosen in the middle of the phylogenetic map which had been obtained from a child with moderate illness to provide an additional comparator within the spectrum of illness severity. The isolates studied included: Methicillin-sensitive S. aureus (MSSA)-29 (mild AHO), Methicillin-resistant S. aureus (MRSA)-12 (moderate AHO), and MRSA-9 (severe AHO). The basis for the selection of these differing clinical strains is to initiate a query into the foundational mechanisms of disease of childhood osteomyelitis. It is our hypothesis that clinical severity of illness is driven by the transcriptional signaling introduced by the pathogenic strains of S.aureus within local communities.

Growth kinetics study and analysis
The clinical isolates were cultured on Tryptic soy agar with 5% sheep blood for 24 hours at 37˚C in ambient air. A pre-inoculum growth culture was prepared by inoculating a single colony into 2 mL of Mueller-Hinton broth (MHB, Thermo Fischer Scientific, Massachusetts, United States; 37˚C, with agitation at 225 rpm, and room air humidity) to achieve a target optical density at 600 nm (OD 600 ) of 0.5. Based on protocol optimization, 0.5 OD 600 occurred at about 4 hours for MSSA-29 (mild), 5 hours for MRSA-12 (moderate), and 5.5 hours for MRSA-9 (severe). For growth kinetics and RNA sequencing culture preparation, an aliquot of 0.5 OD 600 pre-inoculum growth culture was inoculated into 25 mL MHB to a starting OD 600 of 0.025. Cultures were incubated at 37˚C with agitation at 225 rpm and room air humidity. Growth kinetics studies were evaluated by serially removing broth culture aliquots and measuring growth by OD 600 and CFU/mL every hour, for three separate colonies per strain. A one-way ANOVA, followed by the Tukey method for multiple comparisons, was performed to identify statistically significant differences (p<0.05) of OD 600 measurements and rates of change for the three isolates throughout the growth cycle.

RNA isolation, library construction and sequencing
Bacterial RNA isolation was performed from 1.0 mL aliquots of broth culture obtained at 40-minute intervals from 4.5 to 7.8 hours. This time-period correlated with the growth phase intended for transcriptional comparison, with each clinical isolate analyzed in triplicate. RNAstabilized samples were homogenized in 1 ml of MHB and centrifuged at 13,300 x g for 3 minutes at 4˚C to collect the bacterial pellet. RNA was extracted using TRizol reagent (Sigma-Aldrich, United States). The bacterial pellets were then re-suspended in pre-loaded RINO tubes for homogenization of bacterial samples in BBY24M Bullet Blenders1 (Next Advance, United States). Total RNA was extracted using the Bacterial Total RNA extraction Kit (Sigma-Aldrich, USA). Quantity and quality of RNA samples were measured by an Agilent Bioanalyzer 2100 (Agilent Technologies, Santa Clara, California, USA). RNA integrity number (RIN) � 9.0 was considered as quality pass for sequencing. RNA sequencing libraries were prepared with Zymo-Seq RiboFree Total RNA Library Kit (R3000) according to the manufacturer's protocol. Libraries were validated on an Agilent Bioanalyzer 2100. Indexed libraries were equimolarly pooled and sequenced on a SE75 (single-end 75 base pair) Illumina NextSeq550 flow cell. About 20 million sequencing reads were generated for each sample. All sequence data were submitted to NCBI Gene Expression Omnibus (GEO) under the accession number 981186 which is available at the following URL: ID 981186-BioProject-NCBI (nih.gov).

RNA sequence analysis, GO and KEGG pathway enrichment analysis
The NASA RNA-Sequencing Consensus Pipeline (RCP) was used for the analysis [9]. Briefly, the quality of sequence reads was assessed using FastQC/MultiQC. Sequencing reads were mapped to the USA 300 S. aureus reference strain genome of UTSW55 (NCBI GenBank Reference NZ_CP013231.1) [10] using STAR aligner. Counts were then quantified using RSEM. Following this, normalization of counts and identification of differentially expressed genes (DEGs) were performed using DESeq2 [10]. Significant DEGs were selected using fold-change >2 and adjusted p-value <0.05.
Global transcriptional differences between isolates over the progression of growth were evaluated with Principal Component Analysis (PCA) and volcano plots which were generated for three contrasts: MRSA-9 vs MSSA-29 (severe vs mild), MRSA-12 vs MSSA-29 (moderate vs mild), and MRSA-9 vs. MRSA-12 (severe vs moderate). Enrichment of Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways was performed using the Gene Set Enrichment Analysis technique [11,12] with the goal to identify and evaluate similarities and differences of perturbed pathways among the various virulence levels and timepoints. A P-value of < 0.05, adjusted by false discovery rate (FDR), was considered to have statistical significance to achieve significant enrichment.

Bacterial growth kinetics
The three study isolates demonstrated differing growth kinetics based on OD 600 (Fig 1 and Table 1). The mild strain (MSSA-29) had significantly higher OD 600 measurements (p<0.001) at all time points than that of either the moderate (MRSA-12) or severe (MRSA-9) strain. The moderate strain demonstrated a trend of higher OD 600 at 10 hours and significantly greater (p = 0.0013) than the severe isolate at 11 hours of growth. The mild strain reached its peak growth rate in 5 to 6 hours after inoculation, the moderate strain in 8 to 9 hours and the severe strain in 7 to 8 hours.

Cellular aggregation during growth
Polysaccharide intercellular adhesin [13] locus (icaABCD), fibronectin-binding protein (fnbA), fibrinogen-binding protein (fib), and Protein A (spa) comparisons demonstrated significantly differential gene expression between isolates, particularly between the moderate or severe strains when compared with the mild strain ( Table 2). Upregulation of fib, fnb, and ica along with downregulation of spa occurred over the first four time points. The severe and moderate strains showed upregulation of fibrinogen binding of MRSA-9 until the last two time periods of growth when fib, fnbA, and ica became downregulated while spa activity increased. This late (T8-T9) transcriptional pattern was similarly found between MRSA-9 and MSSA-29.  A heat map of the top 50 up-or down-regulated genes (Fig 3) provides greater detail of differences between transcriptional signatures of the isolates. The leukocidins, hemolysins, adhesins and iron metabolism genes are down-regulated for the mild strain from 4.5 to 6.5 hours. These same genes become down-regulated for the severe strain at 7.2 and 7.8 hours. Whereas the genes for lactose metabolism (lac) are highly upregulated in the mild strain during the first four time points and upregulated in the severe strain during the final two time points.

GO and KEGG pathway enrichment analysis
GO pathway analysis illustrates functional differences between isolates. We observed downregulation of biosynthetic, glycolytic, and metabolic processes demonstrated in severe strain as

PLOS ONE
Distinct S. aureus transcriptomes from mild, moderate, and severe osteomyelitis compared to the mild strain. Conversely, upregulated DNA strand exchange and DNA integration was observed over the first four time points of growth in the severe strain (Fig 4). The GO pathway analysis of the moderate versus mild isolates demonstrated similar pattern of transcription to that of the severe versus mild comparison (S1 Fig). However, when the KEGG pathway enrichment analysis showed significant upregulation of SA infection, bacterial invasion, longevity regulating pathways, cell cycle, two component systems, alanine, aspartate and glutamate metabolism, and ABC transporters of the severe strain relative to the mild strain (Fig 5). Whereas functions of metabolism and biosynthesis, including glucose, fructose, galactose, and glycolysis were downregulated in the severe strain. In a similar manner to what was depicted with GO analysis, these changes occurred during the last two periods of growth. KEGG pathway enrichment analysis of moderate to mild and moderate to severe strains are provided separately in supplemental materials (S3 and S4 Figs). Of these comparisons, the differences between severe and moderate strains are minimal during the first four time points. During the last two time points, the severe strain demonstrates upregulated

Discussion
The clinical spectrum of S. aureus osteomyelitis in growing children ranges in severity from mild, with rapid resolution, to severe, with potential for serious long-term adverse outcomes [14]. Although the S. aureus genome is comprised of at least 200 virulence genes, the pathogenetic mechanisms which potentiate virulence and lead to clinical phenotype differentiation, are poorly understood [15][16][17][18]. Panton-Valentine Leucocidin (PVL) presence among S. aureus strains in osteomyelitis has been associated with longer duration of infection and greater inflammatory response [19]. However, it has been difficult to establish a definite connection

PLOS ONE
Distinct S. aureus transcriptomes from mild, moderate, and severe osteomyelitis between a specific virulence gene, such as PVL and clinical disease manifestations given the redundancy of two-component systems and downstream toxins purportedly responsible for the virulence behavior of S. aureus [20]. Previous work has demonstrated genomic heterogeneity among clinical isolates obtained from 71 children with S. aureus AHO which discovered an association between the phylogenetic interconnection of the pathogens and the severity of illness of the affected children. Isolates associated with mild and severe illness were on opposite ends of the phylogenetic map, whereas isolates retrieved from moderately sick children gathered near the center. This suggested the likelihood of bacterial genetic influence upon the clinical manifestations of this disease. To our knowledge, this is the first study to explore transcriptional pathways of S. aureus strains obtained from children with well-differentiated clinical disease based on severity of illness. The findings of this study confirm that growth kinetics and transcriptional activity vary among clinical pathogens obtained from children with well-differentiated disease.
With the emergence of high-throughput technologies of next generation sequencing, transcriptomics, and metabolic network analysis, investigators are increasingly able to explore systems of bacterial metabolic, biosynthetic, and virulence behaviors [21][22][23][24]. One group evaluated the pan-genome of S. aureus using sequence data from 64 strains to identify the core genome, shared by all strains of S. aureus, the accessory genome, present in some, but not all species, and the unique genome, specific to individual isolates [21]. They suggested that S. aureus strains responsible for severe infections may be identified based on growth capabilities and the presence of specific virulence genes [21], which is consistent with our findings. Another group of investigators explored DEGs of infected patients with osteomyelitis using GO and KEGG enrichment analysis of whole blood samples from affected patients and healthy controls [4]. They found 209 SA infection-related genes (SARGs) and 377 osteomyelitisrelated genes (OMRGs) which differentiated infected patients from healthy controls [4]. However, there was no attempt to differentiate clinical phenotypes by severity or explore the underlying pathogenetic mechanisms which may have led to the transcriptional differentiation of the host from healthy control. Important findings of present study are the tendencies for an isolate which caused mild infection (MSSA-29) to favor growth, metabolism, and biosynthesis, while the moderate (MRSA-12) and severe (MRSA-9) isolates favored infection, invasion, two-component systems, quorum sensing, and the ABC transporters. This is supported by the observation of differing growth kinetics in which the mild strain reached the log phase sooner than the moderate or severe strains and then displayed a higher rate of growth. The growth rate and transcriptional differences observed in this study suggest that S. aureus tends to make tradeoffs between growth and virulence.
Another finding of our study is that the bacterial transcription evolved longitudinally over time in an apparent response to environmental changes. As growth progressed and the nutritional resources of the MHB became limited, the growth kinetics and transcription of the isolates evolved accordingly. This is supported by the GO and KEGG pathway findings in which the transcriptional patterns become truncated or reversed during the late time points of growth. It is also notable that the growth and metabolism functions of MRSA-9 and MRSA-12 were significantly down-regulated relative to MSSA-29 prior to the final two time points when these same functions become significantly upregulated in the moderate and severe strains. Other investigators have analyzed events occurring during exponential and post-exponential growth in vitro in S. aureus, showing that virulence gene expression occurs during the growth phase whereas secreted proteins are expressed during the post-exponential phase [19,25,26]. The complexities involved with quorum sensing and transition to alternative transcriptional behavior have been previously studied to demonstrate that codY acts to inhibit metabolic genes with feedback mechanisms triggered by an excess of isoleucine which resulted in growth inhibition in wild type strains but not in codY mutants [27]. In our study, codY was downregulated from T7 to T9 for MRSA-9 relative to MSSA-29.
In the KEGG enrichment pathway analysis the most significantly upregulated function of the severe isolate was that of ABC transporters. This family of membrane proteins participates in ATP-driven transport mechanisms [28]. They hydrolyze ATP, driving a conformational change of transmembrane domains and allow for inward or outward transport of substrates [29,30]. Several ABC transporters have been associated with virulence behavior of S. aureus, including the Ecs ATP-binding cassette shown to be essential for SA expression of the virulence regulatory protein Rot [31,32]. Other ABC transporters implicated in S. aureus pathogenesis include nik, a nickel transporter and cnt, a nickel/cobalt transporter [31,32]. Our data on ABC transporter gene expression in severe strain is consistent with this literature. Another KEGG pathway difference between the moderate and mild strains was that of Staphylococcus aureus infection. S. aureus is known to express several virulence pathways to evade the host immune system, liberate free iron, adhere to bone or collagen and establish infection [31,33,34]. One virulence mechanism includes the modulation of cationic antimicrobial peptides by increasing the positive charge of the cytoplasmic membrane. The two-component system graRS couples an efflux pump (vraG) to an ATPase (vraF) that engages in sensing human defensins [35]. Once these positively charged defensins are detected, S. aureus induces an increase in the surface positive charge and decreases the effective killing of those defensins [35]. Other genes involved in this functional pathway include those encoding for surface proteins like clumping factor B (clfB), iron-regulated surface determinant (isdB) promote adhesion in vitro [36,37]. Another surface protein includes the S. aureus surface protein (sasG) which is known to mask SA microbial surface components recognizing the adhesive matrix molecules (MSCRAMMs) binding to their ligand and promoting biofilm formation [35,38].
Limitations of this study include the inability to extrapolate results to the human condition of childhood osteomyelitis because of the in vitro methodology. The transcriptional differences observed occurred in rich media under standardized growth conditions. Therefore, it is not possible to predict how these isolates would behave in constrained media, an in vivo model, or a human host. However, the results of this study do confirm that organism-specific transcriptional behavior exists among well-differentiated S. aureus clinical isolates. It is therefore plausible that these organism-specific, unique transcriptional signatures would play a role in the clinical phenotype differentiation of children with AHO.
Another study limitation is the comparison of isolates which are growing at different rates through the exponential phase. The decision to isolate RNA longitudinally over time rather than in a manner targeted to specific OD 600 values was intentional. It was discovered during growth curve optimization how challenging it would be to match optical densities. Performing a comparative analysis between isolates at one specific OD was beyond the scope of this study. The data derived from optical density or CFU/mL values does not necessarily reflect the state of active organism growth, given that aggregation may vary substantially between isolates, as was shown in this study. This is additionally confirmed by the finding that the OD 600 of the isolates had only one intersection, approximately 0.75 at T4(4.5h) for MSSA-29, T7(6.5h) for MRSA-12 and T8(7.2h) for MRSA-9. At that intersection of optical densities, there were substantial differences of CFU/mL values with 1.79E+10 for MSSA-29, 3.78E+10 for MRSA-12, and 5.11E+10 for MRSA-9. Investigators have previously shown that S. aureus has a tendency for aggregate formation to occur early in the growth cycle with 50% of cells assembled into aggregates at OD 600 of 0.5 [38]. Key contributors to this tendency to aggregate and form biofilm include S. aureus cell surface proteins (fib and fnb), spa, and PIA (icaABCD) [13,39]. The differential expression of these genes identified in this study lend evidence that aggregation potentially alters growth kinetics. It is notable that the growth of MRSA-12 and MRSA-9 did not differ according to CFU/mL or OD 600 over the six time points of the study. Despite similar growth behavior, these isolates demonstrated unique transcriptional behavior, including that of PIA, fnb and fib during the later stages of growth. Our preliminary growth studies found that MRSA-9 and MRSA-12 significantly diverged at 10 to 11 hours of growth, with the moderate isolate outpacing the severe. The transcriptional distinctions between MRSA-9 and MRSA-12 which occurred at T8 and T9 in this study are premonitory of the kinetic differences that were observed over 11 hours of sustained growth.
While there are numerous DEGs identified in this study, the relative importance of any given pathway to the development or progression of osteomyelitis cannot be ascertained from this data. There are challenges in applying pathway databases of KEGG and GO to a clinical condition that is potentially far more complex that the theoretical network models may be able to depict. Enrichment analysis methods potentially make unrealistic assumptions of statistical independence among genes. Another drawback of pathway enrichment analysis is that these methods ignore genes with no pathway annotations. Ultimately these genes may need to be studied separately in future work.
This study confirms distinct transcriptional profiles among clinical strains of S. aureus in vitro and invokes consideration that the varied transcriptomes may underly the pathogenetic mechanisms leading to well-differentiated phenotypes of disease among children with AHO. Pathogens leading to severe illness may do so by favoring pathways of virulence at the expense of growth being mediated by two-component systems, S. aureus infection, quorum sensing, and ABC transporters. By gaining an understanding of the promotors and inhibitors of these pathways it would be possible to conceive of novel therapeutic strategies to interfere with the mechanisms leading to severe illness. Further work is needed to validate the results of this study using an in vivo model or samples of infected tissue procured during active infections of children with AHO.